A gene’s expression peaks when its first derivative equals zero. Hence, we could use this to order the TFs which have a significant expression peak.
Since the goal of this analysis is ordering of genes, we will ignore multiple testing correction across genes. However, since we want to restrict ourselves to genes that actually do have a significant expression peak, we will do within-gene FWER correction.
Lineage-specific ordering
# neuronal lineage
pdf("../plots/peakAnalysis_Thresholded/neuronal/orderedPeaksPlots.pdf")
neurRes <- orderPeakGenesThresh(grid=grid1, derivatives=derAllTF, tf=tf, lineage=1, sce=sce,
cl=clDatta, main="Neuronal", threshold = 0.1, show_rownames = FALSE)
neurPeaks <- neurRes$firstPeak

dev.off()
## pdf
## 3
pdf("../plots/peakAnalysis_Thresholded/neuronal/orderedPeaksPlots_tallHeatmaps.pdf", width=10, height=65)
neurRes <- orderPeakGenesThresh(grid=grid1, derivatives=derAllTF, tf=tf, lineage=1, sce=sce,
cl=clDatta, main="Neuronal", plotHist=FALSE, plotLines=FALSE,
threshold = 0.1)
dev.off()
## pdf
## 3
# Sus lineage
pdf("../plots/peakAnalysis_Thresholded/sus/orderedPeaksPlots.pdf")
susRes <- orderPeakGenesThresh(grid1=grid2, derivatives=derAllTF, tf=tf, lineage=2, sce=sce,
cl=clDatta, main="Sustentacular", threshold = 0.1)
susPeaks <- susRes$firstPeak
dev.off()
## pdf
## 3
pdf("../plots/peakAnalysis_Thresholded/sus/orderedPeaksPlots_tallHeatmaps.pdf", width=10, height=40)
susRes <- orderPeakGenesThresh(grid1=grid2, derivatives=derAllTF, tf=tf, lineage=2, sce=sce,
cl=clDatta, main="Sustentacular", plotHist=FALSE, plotLines=FALSE, threshold = 0.1)
dev.off()
## pdf
## 3
# rHBC lineage
pdf("../plots/peakAnalysis_Thresholded/sus/orderedPeaksPlots.pdf")
rhbcRes <- orderPeakGenesThresh(grid1=grid3, derivatives=derAllTF, tf=tf, lineage=3, sce=sce,
cl=clDatta, main="rHBC", threshold = 0.1)
rhbcPeaks <- rhbcRes$firstPeak
dev.off()
## pdf
## 3
pdf("../plots/peakAnalysis_Thresholded/rhbc/orderedPeaksPlots_tallHeatmaps.pdf", width=10, height=60)
rhbcRes <- orderPeakGenesThresh(grid1=grid3, derivatives=derAllTF, tf=tf, lineage=3, sce=sce,
cl=clDatta, main="rHBC", plotHist=FALSE, plotLines=FALSE,
threshold = 0.1)
dev.off()
## pdf
## 3
peakList <- list(neur=neurPeaks,
sus=susPeaks,
rhbc=rhbcPeaks)
peakResList <- list(neur=neurRes,
sus=susRes,
rhbc=rhbcRes)
lapply(peakList, length)
## $neur
## [1] 524
##
## $sus
## [1] 231
##
## $rhbc
## [1] 284
saveRDS(peakList, file="../data/peakList_thresholded.rds")
saveRDS(peakResList, file="../data/peakRes_thresholded.rds")
Custom plots for paper
# Neuronal
png("../plots/peakAnalysis_Thresholded/neuronal/lines_neuronal.png", width=9, height=10, units="in", res=300)
neurRes <- orderPeakGenesThresh(grid=grid1, derivatives=derAllTF, tf=tf, lineage=1, sce=sce,
cl=clDatta, main="Neuronal", plotHist=FALSE, plotLines=TRUE,
plotHeatmap = FALSE, threshold = 0.1, clusteredHeatmap = FALSE,
show_rownames=FALSE)
dev.off()
## quartz_off_screen
## 2
png("../plots/peakAnalysis_Thresholded/neuronal/cascadeHeatmap_neuronal.png", width=9, height=10, units="in", res=300)
neurRes <- orderPeakGenesThresh(grid=grid1, derivatives=derAllTF, tf=tf, lineage=1, sce=sce,
cl=clDatta, main="Neuronal", plotHist=FALSE, plotLines=FALSE,
plotHeatmap = TRUE, threshold = 0.1, clusteredHeatmap = FALSE,
show_rownames=FALSE)
dev.off()
## pdf
## 3
# Sus
png("../plots/peakAnalysis_Thresholded/sus/lines_sus.png", width=9, height=10, units="in", res=300)
susRes <- orderPeakGenesThresh(grid=grid2, derivatives=derAllTF, tf=tf, lineage=2, sce=sce,
cl=clDatta, main="Sustentacular", plotHist=FALSE, plotLines=TRUE,
plotHeatmap = FALSE, threshold = 0.1, clusteredHeatmap = FALSE,
show_rownames=FALSE)
dev.off()
## pdf
## 3
png("../plots/peakAnalysis_Thresholded/sus/cascadeHeatmap_sus.png", width=9, height=10, units="in", res=300)
susRes <- orderPeakGenesThresh(grid=grid2, derivatives=derAllTF, tf=tf, lineage=2, sce=sce,
cl=clDatta, main="Sustentacular", plotHist=FALSE, plotLines=FALSE,
plotHeatmap = TRUE, threshold = 0.1, clusteredHeatmap = FALSE,
show_rownames=FALSE)
dev.off()
## pdf
## 3
# HBC
png("../plots/peakAnalysis_Thresholded/rhbc/lines_rhbc.png", width=9, height=10, units="in", res=300)
rhbcRes <- orderPeakGenesThresh(grid=grid3, derivatives=derAllTF, tf=tf, lineage=3, sce=sce,
cl=clDatta, main="rHBC", plotHist=FALSE, plotLines=TRUE,
plotHeatmap = FALSE, threshold = 0.1, clusteredHeatmap = FALSE,
show_rownames=FALSE)
dev.off()
## pdf
## 3
png("../plots/peakAnalysis_Thresholded/rhbc/cascadeHeatmap_rhbc.png", width=9, height=10, units="in", res=300)
rhbcRes <- orderPeakGenesThresh(grid=grid3, derivatives=derAllTF, tf=tf, lineage=3, sce=sce,
cl=clDatta, main="rHBC", plotHist=FALSE, plotLines=FALSE,
plotHeatmap = TRUE, threshold = 0.1, clusteredHeatmap = FALSE,
show_rownames=FALSE)
dev.off()
## pdf
## 3
All cascades next to each other
phNeur <- orderPeakGenesThresh(grid=grid1, derivatives=derAllTF, tf=tf, lineage=1, sce=sce,
cl=clDatta, main="Neuronal", plotHist=FALSE, plotLines=FALSE,
plotHeatmap = TRUE, threshold = 0.1, clusteredHeatmap = FALSE,
show_rownames=FALSE, returnOnlyHeatmap = TRUE)

phSus <- orderPeakGenesThresh(grid=grid2, derivatives=derAllTF, tf=tf, lineage=2, sce=sce,
cl=clDatta, main="Sustentacular", plotHist=FALSE, plotLines=FALSE,
plotHeatmap = TRUE, threshold = 0.1, clusteredHeatmap = FALSE,
show_rownames=FALSE, returnOnlyHeatmap = TRUE)

phHBC <- orderPeakGenesThresh(grid=grid3, derivatives=derAllTF, tf=tf, lineage=3, sce=sce,
cl=clDatta, main="rHBC", plotHist=FALSE, plotLines=FALSE,
plotHeatmap = TRUE, threshold = 0.1, clusteredHeatmap = FALSE,
show_rownames=FALSE, returnOnlyHeatmap = TRUE)

cowplot::plot_grid(phNeur[[4]], phSus[[4]], phHBC[[4]],
nrow=1, ncol=3)

ggsave("../plots/fig2_allCascades.png")
## Saving 7 x 5 in image
Identify shared and unique lineage TFs
We can use the results above to identify shared and unique TFs for each lineage, using a set diff. For shared TFs, we first aggregate the p-values testing whether a derivative is different from zero within the gene, using Sidak aggregation. Next, we perform FDR correction on the aggregated p-values for each lineage separately. The TFs that are significant post FDR correction in all three lineages are considered to be shared.
For TFs uniquely peaking in a lineage, we retain TFs that are significantly peaking in that lineage, and that are higher expressed in that lineage versus the other two lineages, by a fold change threshold of \(1.5\) on the maximum fitted expression value.
yhatAllTF <- predictSmooth(sce, gene=tf, nPoints=100, tidy=FALSE)
library(aggregation)
## shared TFs
# aggregate p-values using Sidak method for each lineage
aggPvals <- cbind(apply(neurRes$pvalThresh, 1, sidak),
apply(rhbcRes$pvalThresh, 1, sidak),
apply(susRes$pvalThresh, 1, sidak))
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
## Warning in FUN(newX[, i], ...): Extreme p-values around and below 10e-17 will
## produce aggregated p-value of 0. Replace extreme p-values with 10e-17 to obtain
## upper bound on the aggregated p-value
# FDR correction for each lineage across all TFs
aggPvalsFDR <- apply(aggPvals,2,p.adjust,method="fdr")
barplot(table(rowSums(aggPvalsFDR <= 0.05)), xlab="Number of lineages in which a TF significantly peaks.")

sharedTF <- rownames(aggPvalsFDR)[which(rowSums(aggPvalsFDR <= 0.05) == 3)]
length(sharedTF)
## [1] 65
## unique TFs
getUniqueTFs <- function(peakResList, yhatAllTF, lineage, name, fcThreshold=1.5){
# get TFs identified for that lineage
tfLineage <- names(peakResList[[name]]$firstPeak)
# within that set, check if expression is much higher
linID <- ((lineage-1)*100+1):(lineage*100)
delta <- rowMax(yhatAllTF[tfLineage,linID]) / rowMax(yhatAllTF[tfLineage,-linID])
tfLineageDelta <- tfLineage[delta > fcThreshold]
return(tfLineageDelta)
}
neurUniqTFs <- getUniqueTFs(peakResList, yhatAllTF, lineage=1, name="neur")
susUniqTFs <- getUniqueTFs(peakResList, yhatAllTF, lineage=2, name="sus")
rhbcUniqTFs <- getUniqueTFs(peakResList, yhatAllTF, lineage=3, name="rhbc")
length(neurUniqTFs) ; length(susUniqTFs) ; length(rhbcUniqTFs)
## [1] 352
## [1] 31
## [1] 61
# write.table(sharedTF, file="../data/sharedTF.txt", row.names=FALSE, col.names=FALSE, quote=FALSE)
# write.table(neurUniqTFs, file="../data/neurUniqTf_threshold.txt", row.names=FALSE, col.names=FALSE, quote=FALSE)
# write.table(susUniqTFs, file="../data/susUniqTf_threshold.txt", row.names=FALSE, col.names=FALSE, quote=FALSE)
# write.table(rhbcUniqTFs, file="../data/rhbcUniqTf_threshold.txt", row.names=FALSE, col.names=FALSE, quote=FALSE)
Alternative approach to identify shared TFs using Fisher p-value aggregation
yhatAllTF <- predictSmooth(sce, gene=tf, nPoints=100, tidy=FALSE)
fisherX <- function(pval){
pval[pval==0] <- 1e-100
stat <- -2 * sum(log(pval))
return(stat)
}
## shared TFs
# aggregate p-values using Fisher method for each lineage
aggStats <- cbind(apply(neurRes$pvalThresh, 1, fisherX),
apply(rhbcRes$pvalThresh, 1, fisherX),
apply(susRes$pvalThresh, 1, fisherX))
plot(density(log(aggStats[,1])), ylim=c(0,0.4))
lines(density(log(aggStats[,2])), col="blue")
lines(density(log(aggStats[,3])), col="darkseagreen3")
abline(v=7.5, col="red")

sharedTF1 <- rownames(aggStats)[which(rowSums(log(aggStats) > 7.5) == 3)]
sharedTF1
## [1] "Ebf1" "Lhx2" "Fos" "Nupr1"
## unique TFs
getUniqueTFs <- function(peakResList, yhatAllTF, lineage, name, fcThreshold=1.5){
# get TFs identified for that lineage
tfLineage <- names(peakResList[[name]]$firstPeak)
# within that set, check if expression is much higher
linID <- ((lineage-1)*100+1):(lineage*100)
delta <- rowMax(yhatAllTF[tfLineage,linID]) / rowMax(yhatAllTF[tfLineage,-linID])
tfLineageDelta <- tfLineage[delta > fcThreshold]
return(tfLineageDelta)
}
neurUniqTFs <- getUniqueTFs(peakResList, yhatAllTF, lineage=1, name="neur")
susUniqTFs <- getUniqueTFs(peakResList, yhatAllTF, lineage=2, name="sus")
rhbcUniqTFs <- getUniqueTFs(peakResList, yhatAllTF, lineage=3, name="rhbc")
length(neurUniqTFs) ; length(susUniqTFs) ; length(rhbcUniqTFs)
## [1] 352
## [1] 31
## [1] 61
write.table(sharedTF1, file="../data/sharedTF_v2.txt", row.names=FALSE, col.names=FALSE, quote=FALSE)
write.table(neurUniqTFs, file="../data/neurUniqTf_threshold.txt", row.names=FALSE, col.names=FALSE, quote=FALSE)
write.table(susUniqTFs, file="../data/susUniqTf_threshold.txt", row.names=FALSE, col.names=FALSE, quote=FALSE)
write.table(rhbcUniqTFs, file="../data/rhbcUniqTf_threshold.txt", row.names=FALSE, col.names=FALSE, quote=FALSE)
Another alternative approach, using pseudotime agreement
# 90 TFs in all three
sharedTF1 <- intersect(intersect(names(neurRes$firstPeak), names(susRes$firstPeak)), names(rhbcRes$firstPeak))
# check for similar pseudotime
ptSharedTF1 <- data.frame(neur=neurRes$firstPeak[sharedTF1],
sus=susRes$firstPeak[sharedTF1],
rhbc=rhbcRes$firstPeak[sharedTF1])
deltaTF1 <- rowDiffs(rowRanges(as.matrix(ptSharedTF1)))
hist(deltaTF1, breaks=50)

plot(density(log(deltaTF1)))

sharedTFDelta <- sharedTF1[deltaTF1[,1] < 1]
heatmapScaledAcrossAllLineages(models=sce,
genes=sharedTFDelta,
nPoints=100,
sds=sds,
cl=clDatta,
height=5,
width=5,
showRowNames = FALSE,
showLegend = FALSE,
outFile="../plots/peakAnalysis_Thresholded/shared/sharedTF_scaledAcrossAllLineages_delta.pdf")



heatmapScaledAcrossAllLineages(models=sce,
genes=sharedTFDelta,
nPoints=100,
sds=sds,
cl=clDatta,
height=7,
width=5,
showRowNames = TRUE,
showLegend = FALSE,
outFile="../plots/peakAnalysis_Thresholded/shared/sharedTF_scaledAcrossAllLineages_delta_tfNames.pdf")


